##Examen des Données:
Avant de se plonger dans des analyses spécifiques, il est essentiel d’examiner les données pour comprendre leur structure, leurs variables et leurs statistiques de base. Nous chargerons les données, examinerons les premières lignes et passerons en revue les types de données et les statistiques sommaires.
# Affichage des premières lignes de l'ensemble de données journalier
head(day_data)
## # A tibble: 6 × 16
## instant dteday season yr mnth holiday weekday workingday weathersit
## <dbl> <date> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 1 2011-01-01 1 0 1 0 6 0 2
## 2 2 2011-01-02 1 0 1 0 0 0 2
## 3 3 2011-01-03 1 0 1 0 1 1 1
## 4 4 2011-01-04 1 0 1 0 2 1 1
## 5 5 2011-01-05 1 0 1 0 3 1 1
## 6 6 2011-01-06 1 0 1 0 4 1 1
## # ℹ 7 more variables: temp <dbl>, atemp <dbl>, hum <dbl>, windspeed <dbl>,
## # casual <dbl>, registered <dbl>, cnt <dbl>
# Résumé de l'ensemble de données journalier
summary(day_data)
## instant dteday season yr
## Min. : 1.0 Min. :2011-01-01 Min. :1.000 Min. :0.0000
## 1st Qu.:183.5 1st Qu.:2011-07-02 1st Qu.:2.000 1st Qu.:0.0000
## Median :366.0 Median :2012-01-01 Median :3.000 Median :1.0000
## Mean :366.0 Mean :2012-01-01 Mean :2.497 Mean :0.5007
## 3rd Qu.:548.5 3rd Qu.:2012-07-01 3rd Qu.:3.000 3rd Qu.:1.0000
## Max. :731.0 Max. :2012-12-31 Max. :4.000 Max. :1.0000
## mnth holiday weekday workingday
## Min. : 1.00 Min. :0.00000 Min. :0.000 Min. :0.000
## 1st Qu.: 4.00 1st Qu.:0.00000 1st Qu.:1.000 1st Qu.:0.000
## Median : 7.00 Median :0.00000 Median :3.000 Median :1.000
## Mean : 6.52 Mean :0.02873 Mean :2.997 Mean :0.684
## 3rd Qu.:10.00 3rd Qu.:0.00000 3rd Qu.:5.000 3rd Qu.:1.000
## Max. :12.00 Max. :1.00000 Max. :6.000 Max. :1.000
## weathersit temp atemp hum
## Min. :1.000 Min. :0.05913 Min. :0.07907 Min. :0.0000
## 1st Qu.:1.000 1st Qu.:0.33708 1st Qu.:0.33784 1st Qu.:0.5200
## Median :1.000 Median :0.49833 Median :0.48673 Median :0.6267
## Mean :1.395 Mean :0.49538 Mean :0.47435 Mean :0.6279
## 3rd Qu.:2.000 3rd Qu.:0.65542 3rd Qu.:0.60860 3rd Qu.:0.7302
## Max. :3.000 Max. :0.86167 Max. :0.84090 Max. :0.9725
## windspeed casual registered cnt
## Min. :0.02239 Min. : 2.0 Min. : 20 Min. : 22
## 1st Qu.:0.13495 1st Qu.: 315.5 1st Qu.:2497 1st Qu.:3152
## Median :0.18097 Median : 713.0 Median :3662 Median :4548
## Mean :0.19049 Mean : 848.2 Mean :3656 Mean :4504
## 3rd Qu.:0.23321 3rd Qu.:1096.0 3rd Qu.:4776 3rd Qu.:5956
## Max. :0.50746 Max. :3410.0 Max. :6946 Max. :8714
# Structure de l'ensemble de données journalier
str(day_data)
## spc_tbl_ [731 × 16] (S3: spec_tbl_df/tbl_df/tbl/data.frame)
## $ instant : num [1:731] 1 2 3 4 5 6 7 8 9 10 ...
## $ dteday : Date[1:731], format: "2011-01-01" "2011-01-02" ...
## $ season : num [1:731] 1 1 1 1 1 1 1 1 1 1 ...
## $ yr : num [1:731] 0 0 0 0 0 0 0 0 0 0 ...
## $ mnth : num [1:731] 1 1 1 1 1 1 1 1 1 1 ...
## $ holiday : num [1:731] 0 0 0 0 0 0 0 0 0 0 ...
## $ weekday : num [1:731] 6 0 1 2 3 4 5 6 0 1 ...
## $ workingday: num [1:731] 0 0 1 1 1 1 1 0 0 1 ...
## $ weathersit: num [1:731] 2 2 1 1 1 1 2 2 1 1 ...
## $ temp : num [1:731] 0.344 0.363 0.196 0.2 0.227 ...
## $ atemp : num [1:731] 0.364 0.354 0.189 0.212 0.229 ...
## $ hum : num [1:731] 0.806 0.696 0.437 0.59 0.437 ...
## $ windspeed : num [1:731] 0.16 0.249 0.248 0.16 0.187 ...
## $ casual : num [1:731] 331 131 120 108 82 88 148 68 54 41 ...
## $ registered: num [1:731] 654 670 1229 1454 1518 ...
## $ cnt : num [1:731] 985 801 1349 1562 1600 ...
## - attr(*, "spec")=
## .. cols(
## .. instant = col_double(),
## .. dteday = col_date(format = ""),
## .. season = col_double(),
## .. yr = col_double(),
## .. mnth = col_double(),
## .. holiday = col_double(),
## .. weekday = col_double(),
## .. workingday = col_double(),
## .. weathersit = col_double(),
## .. temp = col_double(),
## .. atemp = col_double(),
## .. hum = col_double(),
## .. windspeed = col_double(),
## .. casual = col_double(),
## .. registered = col_double(),
## .. cnt = col_double()
## .. )
## - attr(*, "problems")=<externalptr>
## # A tibble: 6 × 17
## instant dteday season yr mnth hr holiday weekday workingday
## <dbl> <date> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 1 2011-01-01 1 0 1 0 0 6 0
## 2 2 2011-01-01 1 0 1 1 0 6 0
## 3 3 2011-01-01 1 0 1 2 0 6 0
## 4 4 2011-01-01 1 0 1 3 0 6 0
## 5 5 2011-01-01 1 0 1 4 0 6 0
## 6 6 2011-01-01 1 0 1 5 0 6 0
## # ℹ 8 more variables: weathersit <dbl>, temp <dbl>, atemp <dbl>, hum <dbl>,
## # windspeed <dbl>, casual <dbl>, registered <dbl>, cnt <dbl>
## instant dteday season yr
## Min. : 1 Min. :2011-01-01 Min. :1.000 Min. :0.0000
## 1st Qu.: 4346 1st Qu.:2011-07-04 1st Qu.:2.000 1st Qu.:0.0000
## Median : 8690 Median :2012-01-02 Median :3.000 Median :1.0000
## Mean : 8690 Mean :2012-01-02 Mean :2.502 Mean :0.5026
## 3rd Qu.:13034 3rd Qu.:2012-07-02 3rd Qu.:3.000 3rd Qu.:1.0000
## Max. :17379 Max. :2012-12-31 Max. :4.000 Max. :1.0000
## mnth hr holiday weekday
## Min. : 1.000 Min. : 0.00 Min. :0.00000 Min. :0.000
## 1st Qu.: 4.000 1st Qu.: 6.00 1st Qu.:0.00000 1st Qu.:1.000
## Median : 7.000 Median :12.00 Median :0.00000 Median :3.000
## Mean : 6.538 Mean :11.55 Mean :0.02877 Mean :3.004
## 3rd Qu.:10.000 3rd Qu.:18.00 3rd Qu.:0.00000 3rd Qu.:5.000
## Max. :12.000 Max. :23.00 Max. :1.00000 Max. :6.000
## workingday weathersit temp atemp
## Min. :0.0000 Min. :1.000 Min. :0.020 Min. :0.0000
## 1st Qu.:0.0000 1st Qu.:1.000 1st Qu.:0.340 1st Qu.:0.3333
## Median :1.0000 Median :1.000 Median :0.500 Median :0.4848
## Mean :0.6827 Mean :1.425 Mean :0.497 Mean :0.4758
## 3rd Qu.:1.0000 3rd Qu.:2.000 3rd Qu.:0.660 3rd Qu.:0.6212
## Max. :1.0000 Max. :4.000 Max. :1.000 Max. :1.0000
## hum windspeed casual registered
## Min. :0.0000 Min. :0.0000 Min. : 0.00 Min. : 0.0
## 1st Qu.:0.4800 1st Qu.:0.1045 1st Qu.: 4.00 1st Qu.: 34.0
## Median :0.6300 Median :0.1940 Median : 17.00 Median :115.0
## Mean :0.6272 Mean :0.1901 Mean : 35.68 Mean :153.8
## 3rd Qu.:0.7800 3rd Qu.:0.2537 3rd Qu.: 48.00 3rd Qu.:220.0
## Max. :1.0000 Max. :0.8507 Max. :367.00 Max. :886.0
## cnt
## Min. : 1.0
## 1st Qu.: 40.0
## Median :142.0
## Mean :189.5
## 3rd Qu.:281.0
## Max. :977.0
## spc_tbl_ [17,379 × 17] (S3: spec_tbl_df/tbl_df/tbl/data.frame)
## $ instant : num [1:17379] 1 2 3 4 5 6 7 8 9 10 ...
## $ dteday : Date[1:17379], format: "2011-01-01" "2011-01-01" ...
## $ season : num [1:17379] 1 1 1 1 1 1 1 1 1 1 ...
## $ yr : num [1:17379] 0 0 0 0 0 0 0 0 0 0 ...
## $ mnth : num [1:17379] 1 1 1 1 1 1 1 1 1 1 ...
## $ hr : num [1:17379] 0 1 2 3 4 5 6 7 8 9 ...
## $ holiday : num [1:17379] 0 0 0 0 0 0 0 0 0 0 ...
## $ weekday : num [1:17379] 6 6 6 6 6 6 6 6 6 6 ...
## $ workingday: num [1:17379] 0 0 0 0 0 0 0 0 0 0 ...
## $ weathersit: num [1:17379] 1 1 1 1 1 2 1 1 1 1 ...
## $ temp : num [1:17379] 0.24 0.22 0.22 0.24 0.24 0.24 0.22 0.2 0.24 0.32 ...
## $ atemp : num [1:17379] 0.288 0.273 0.273 0.288 0.288 ...
## $ hum : num [1:17379] 0.81 0.8 0.8 0.75 0.75 0.75 0.8 0.86 0.75 0.76 ...
## $ windspeed : num [1:17379] 0 0 0 0 0 0.0896 0 0 0 0 ...
## $ casual : num [1:17379] 3 8 5 3 0 0 2 1 1 8 ...
## $ registered: num [1:17379] 13 32 27 10 1 1 0 2 7 6 ...
## $ cnt : num [1:17379] 16 40 32 13 1 1 2 3 8 14 ...
## - attr(*, "spec")=
## .. cols(
## .. instant = col_double(),
## .. dteday = col_date(format = ""),
## .. season = col_double(),
## .. yr = col_double(),
## .. mnth = col_double(),
## .. hr = col_double(),
## .. holiday = col_double(),
## .. weekday = col_double(),
## .. workingday = col_double(),
## .. weathersit = col_double(),
## .. temp = col_double(),
## .. atemp = col_double(),
## .. hum = col_double(),
## .. windspeed = col_double(),
## .. casual = col_double(),
## .. registered = col_double(),
## .. cnt = col_double()
## .. )
## - attr(*, "problems")=<externalptr>
sum(is.na(hour_data))
## [1] 0
sum(is.na(day_data))
## [1] 0
sum(duplicated(hour_data))
## [1] 0
sum(duplicated(day_data))
## [1] 0
Nous allons maintenant examiner nos données pour répondre aux questions spécifiques posées.
##Changements de Température Selon les Saisons #Comment les températures changent-elles selon les saisons ? Quelles sont les températures moyennes et médianes ?
Pour répondre à cette question, nous analyserons la température (temp) par saison et calculerons les moyennes et médianes.
library(dplyr)
# Calcul de la température moyenne et médiane par saison
temperature_stats <- day_data %>%
group_by(season) %>%
summarise(mean_temp = mean(temp, na.rm = TRUE),
median_temp = median(temp, na.rm = TRUE))
# Affichage des résultats
print(temperature_stats)
## # A tibble: 4 × 3
## season mean_temp median_temp
## <dbl> <dbl> <dbl>
## 1 1 0.298 0.286
## 2 2 0.544 0.562
## 3 3 0.706 0.715
## 4 4 0.423 0.409
Cette partie du code permet d’obtenir une vue d’ensemble des températures par saison, ce qui est crucial pour comprendre comment les conditions météorologiques influencent l’utilisation des vélos.
##Y a-t-il une corrélation entre la température (temp/atemp) et le nombre total de locations de vélos ?
Nous examinerons la corrélation entre la température, la température ressentie et le nombre total de locations.
# Charger la bibliothèque corrplot si ce n'est pas déjà fait
# Si vous ne l'avez pas déjà installée, exécutez install.packages("corrplot")
library(corrplot)
## corrplot 0.92 loaded
# Calculer la température moyenne entre temp et atemp
day_data$mean_temp_atemp <- rowMeans(day_data[, c("temp", "atemp")])
# Calculer la matrice de corrélation
correlation_matrix <- cor(day_data[, c("temp", "atemp", "mean_temp_atemp", "cnt")], use = "complete.obs")
# Visualiser la matrice de corrélation
corrplot(correlation_matrix, method = "circle")
#Quelles sont les températures moyennes, l’humidité, la vitesse du vent et les locations totales par mois ?
Nous calculerons ces moyennes pour chaque mois pour identifier des tendances saisonnières.
Les valeurs de température, d’humidité et de vitesse du vent sont normalisées. La température moyenne est la plus élevée en juillet (0.755) et la plus basse en décembre (0.324). L’humidité est relativement stable tout au long de l’année avec une légère augmentation pendant les mois d’été. La vitesse du vent moyenne est généralement la plus faible en juillet (0.166) et la plus élevée en avril (0.234). Le nombre de locations de vélos est le plus élevé en juin (5772.37) et le plus bas en janvier (2176.34).
# Calcul des moyennes mensuelles
monthly_averages <- day_data %>%
group_by(mnth) %>%
summarise(
mean_temp = mean(temp, na.rm = TRUE),
mean_humidity = mean(hum, na.rm = TRUE),
mean_windspeed = mean(windspeed, na.rm = TRUE),
mean_rentals = mean(cnt, na.rm = TRUE)
)
# Affichage des moyennes mensuelles
print(monthly_averages)
## # A tibble: 12 × 5
## mnth mean_temp mean_humidity mean_windspeed mean_rentals
## <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 1 0.236 0.586 0.206 2176.
## 2 2 0.299 0.567 0.216 2655.
## 3 3 0.391 0.588 0.223 3692.
## 4 4 0.470 0.588 0.234 4485.
## 5 5 0.595 0.689 0.183 5350.
## 6 6 0.684 0.576 0.185 5772.
## 7 7 0.755 0.598 0.166 5564.
## 8 8 0.709 0.638 0.173 5664.
## 9 9 0.616 0.715 0.166 5767.
## 10 10 0.485 0.694 0.175 5199.
## 11 11 0.369 0.625 0.184 4247.
## 12 12 0.324 0.666 0.177 3404.
#La température est-elle associée aux locations de vélos (enregistrés vs occasionnels) ?
Nous étudierons comment la température affecte différemment les utilisateurs enregistrés et occasionnels.
# Association de la température avec les locations de vélos
rentals_by_temp <- day_data %>%
group_by(temp) %>%
summarise(mean_casual = mean(casual),
mean_registered = mean(registered))
# Visualisation
ggplot(rentals_by_temp, aes(x = temp)) +
geom_line(aes(y = mean_casual, color = "Occasionnel")) +
geom_line(aes(y = mean_registered, color = "Enregistré")) +
labs(title = "Locations de Vélos par Température",
x = "Température Normalisée",
y = "Locations Moyennes")
# Calcul de la corrélation entre température et locations pour les utilisateurs enregistrés et occasionnels
correlation_temp_users <- day_data %>%
select(temp, casual, registered) %>%
cor(use = "complete.obs")
# Affichage de la matrice de corrélation
print(correlation_temp_users)
## temp casual registered
## temp 1.0000000 0.5432847 0.5400120
## casual 0.5432847 1.0000000 0.3952825
## registered 0.5400120 0.3952825 1.0000000
Calculer la corrélation entre la température et les locations pour les utilisateurs enregistrés et occasionnels.
Résumer comment les locations moyennes varient avec la température en regroupant les températures par intervalles puis en calculant le nombre moyen de locations pour chaque groupe. ????????????????????je ne sais pas
# Conversion de 'dteday' en format de date
day_data$dteday <- as.Date(day_data$dteday, format="%Y-%m-%d")
# Tracé de 'cnt' vs 'dteday'
ggplot(day_data, aes(x=dteday, y=cnt)) +
geom_line() +
labs(title="Nombre de locations de vélos par jour",
x="Date",
y="Nombre total de locations de vélos") +
theme(axis.text.x = element_text(angle=90, hjust=1))
Il y a une saison annuelle, on voit une tendance ( mois plus chauds
implique hausse des locations et une baisse pendant les mois plus
froids). La variance est plus au moins constante.
# Conversion de 'dteday' en format de date
day_data$dteday <- as.Date(day_data$dteday, format="%Y-%m-%d")
# Identifier les valeurs aberrantes, par exemple en utilisant la méthode des écarts interquartiles
IQR_values <- IQR(day_data$cnt)
quantiles <- quantile(day_data$cnt, probs=c(.25, .75))
cap <- quantiles[2] + 1.5 * IQR_values
floor <- quantiles[1] - 1.5 * IQR_values
# Filtrer les valeurs aberrantes
day_data <- day_data %>%
filter(cnt >= floor & cnt <= cap)
# Lisser la série temporelle en utilisant une moyenne mobile avec la fonction 'filter' du package 'stats'
day_data$cnttimeseries <- ts(day_data$cnt, frequency = 12)
hw_model <- HoltWinters(day_data$cnttimeseries)
plot(hw_model)
# Extraire les valeurs prévues à partir du modèle HoltWinters
fitted_values <- hw_model$fitted[, "xhat"]
# Ajouter la fréquence appropriée à la série lissée (remplacez 12 par la fréquence appropriée)
smoothed_ts <- ts(fitted_values, frequency = 12)
plot(smoothed_ts)
result <- adf.test(smoothed_ts, alternative = "stationary")
# Afficher les résultats du test
print(result)
##
## Augmented Dickey-Fuller Test
##
## data: smoothed_ts
## Dickey-Fuller = -0.7527, Lag order = 8, p-value = 0.9658
## alternative hypothesis: stationary
La série n’est pas stationnaire (variance constante mais moyenne variable) le test de Dicker Fuller confirme cela avec une p_value = 0.96 (hyp nulle : non stationnaire), on distingue une saison (qui se répète deux fois).
# Appliquer la différenciation saisonnière
diff_series <- diff(smoothed_ts, lag=12)
# Tracer la série différenciée
plot(diff_series, main="Série Temporelle Différenciée", xlab="Temps", ylab="Différences")
pacf(diff_series,lag.max=100, main="PACF")
acf(diff_series,lag.max= 100, main="ACF")
On remarque le PACF décroit exponentiellement, et sur le ACF on lit :
MA(9) On remarque le ACF décroit exponentiellement, et sur le PACF on
lit : AR(2)
Concernant les saisons, sur le ACF ça décroit et sur le PACF on lit P = 4 sur le PACF ça décroit et sur le ACF on lit Q=4
#MA(9) ou AR(2)
#ARIMA(2,0,0)(3,1,0)h=12 ***
#ARIMA(0,0,9)(4,1,0)h=12 ****
#ARIMA(2,0,0)(0,1,2)h=12 *
#ARIMA(0,0,9)(0,1,2)h=12 **
model_m1 <- arima(smoothed_ts, order = c(2, 0, 0), seasonal = list(order = c(3, 1, 0), period = 12))
residuals_m1 <- residuals(model_m1)
plot(residuals_m1)
#mean=0 var=cte
par(mfrow = c(2, 1))
acf(residuals_m1, main = "ACF of Residuals")
pacf(residuals_m1, main = "PACF of Residuals")
#no significant peak on the ACF and PACF
Box.test(residuals_m1, type = "Ljung-Box")
##
## Box-Ljung test
##
## data: residuals_m1
## X-squared = 0.057693, df = 1, p-value = 0.8102
#p-value= 0.8102
aic_value_m1 <- AIC(model_m1)
cat("AIC du modèle m1:", aic_value_m1, "\n")
## AIC du modèle m1: 10179.9
#AIC=10179.9
model_m2 <- arima(smoothed_ts, order = c(0, 0, 9), seasonal = list(order = c(4, 1, 0), period = 12))
residuals_m2 <- residuals(model_m2)
plot(residuals_m2)
#mean=0 var=cte
par(mfrow = c(2, 1))
acf(residuals_m2, main = "ACF of Residuals")
pacf(residuals_m2, main = "PACF of Residuals")
#no significant peak on the ACF and PACF
Box.test(residuals_m2, type = "Ljung-Box")
##
## Box-Ljung test
##
## data: residuals_m2
## X-squared = 0.064682, df = 1, p-value = 0.7992
#p-value= 0.7992
aic_value_m2 <- AIC(model_m2)
cat("AIC du modèle m2:", aic_value_m2, "\n")
## AIC du modèle m2: 10213.16
#AIC=10213.16
model_m3 <- arima(smoothed_ts, order = c(2, 0, 0), seasonal = list(order = c(0, 1, 2), period = 12))
residuals_m3 <- residuals(model_m3)
plot(residuals_m3)
#mean=0 var=cte
par(mfrow = c(2, 1))
acf(residuals_m3, main = "ACF of Residuals")
pacf(residuals_m3, main = "PACF of Residuals")
#no significant peak on the ACF and PACF
Box.test(residuals_m3, type = "Ljung-Box")
##
## Box-Ljung test
##
## data: residuals_m3
## X-squared = 0.126, df = 1, p-value = 0.7226
#p-value= 0.7226
aic_value_m3 <- AIC(model_m3)
cat("AIC du modèle m3:", aic_value_m3, "\n")
## AIC du modèle m3: 10131.05
#AIC=10131.05
model_m4 <- arima(smoothed_ts, order = c(0, 0, 9), seasonal = list(order = c(0, 1, 2), period = 12))
residuals_m4 <- residuals(model_m4)
plot(residuals_m4)
#mean=0 var=cte
par(mfrow = c(2, 1))
acf(residuals_m4, main = "ACF of Residuals")
pacf(residuals_m4, main = "PACF of Residuals")
#no significant peak on the ACF and PACF
Box.test(residuals_m4, type = "Ljung-Box")
##
## Box-Ljung test
##
## data: residuals_m4
## X-squared = 0.051688, df = 1, p-value = 0.8202
#p-value= 0.8202
aic_value_m4 <- AIC(model_m4)
cat("AIC du modèle m4:", aic_value_m4, "\n")
## AIC du modèle m4: 10206.35
#AIC=10206.35
Le meilleur modèle a l’air d’être le modèle m2 : ARIMA(0,0,9)(4,1,0)h=12
library(forecast)
# Ajuster un modèle ARIMA automatiquement sur les données stationnaires
# Supposons que 'diff_log_cnt_clean' est votre série temporelle prétraitée
model <- auto.arima(smoothed_ts)
# Afficher le résumé du modèle ajusté
summary(model)
## Series: smoothed_ts
## ARIMA(5,1,1)(2,0,0)[12]
##
## Coefficients:
## ar1 ar2 ar3 ar4 ar5 ma1 sar1 sar2
## -0.5396 -0.0782 -0.1907 -0.2164 -0.1897 0.6284 0.2808 0.3158
## s.e. 0.2044 0.0463 0.0472 0.0470 0.0402 0.2054 0.0372 0.0370
##
## sigma^2 = 104406: log likelihood = -5165.91
## AIC=10349.81 AICc=10350.07 BIC=10391
##
## Training set error measures:
## ME RMSE MAE MPE MAPE MASE
## Training set -3.822976 321.0912 235.8936 -0.6894523 6.372798 0.4056763
## ACF1
## Training set 0.002163913
##Désaisonnaliser les données :
#Désaisonnaliser les données
diff_cnt <- diff(day_data$cnt, lag=12)
plot(diff_cnt)
pacf(diff_cnt, main="PACF")
L’ACF présente une décroissance rapide après le premier décalage, ce qui
suggère un terme MA. Si ce décalage est le seul significatif, cela
pourrait indiquer un modèle MA(1).
Le PACF montre un pic significatif au premier décalage et ensuite il se stabilise, ce qui est typique d’un modèle AR(1).
acf(diff_cnt, main="ACF")
#AR(1)
#MA(3)
#ARIMA(1,0,0)
#ARIMA(0,0,3)
model_m1_deseasonal <- arima(diff_cnt, order = c(1, 0, 0))
residuals_m1_deseasonal <- residuals(model_m1_deseasonal)
plot(residuals_m1_deseasonal)
#mean=0 var=cte
par(mfrow = c(2, 1))
acf(residuals_m1_deseasonal, main = "ACF of Residuals")
pacf(residuals_m1_deseasonal, main = "PACF of Residuals")
#no significant peak on the ACF and PACF
Box.test(residuals_m1_deseasonal, type = "Ljung-Box")
##
## Box-Ljung test
##
## data: residuals_m1_deseasonal
## X-squared = 0.28091, df = 1, p-value = 0.5961
#p-value= 0.5961
aic_value_m1_deseasonal<- AIC(model_m1_deseasonal)
cat("AIC du modèle m4:", aic_value_m1_deseasonal, "\n")
## AIC du modèle m4: 12343.73
#AIC=12343.73
model_m2_deseasonal <- arima(diff_cnt, order = c(0, 0, 3))
residuals_m2_deseasonal <- residuals(model_m2_deseasonal)
plot(residuals_m2_deseasonal)
#mean=0 var=cte
par(mfrow = c(2, 1))
acf(residuals_m2_deseasonal, main = "ACF of Residuals")
pacf(residuals_m2_deseasonal, main = "PACF of Residuals")
#no significant peak on the ACF and PACF
Box.test(residuals_m2_deseasonal, type = "Ljung-Box")
##
## Box-Ljung test
##
## data: residuals_m2_deseasonal
## X-squared = 0.0016148, df = 1, p-value = 0.9679
#p-value= 0.9679
aic_value_m2_deseasonal<- AIC(model_m2_deseasonal)
cat("AIC du modèle m2:", aic_value_m2_deseasonal, "\n")
## AIC du modèle m2: 12346.61
#AIC=12346.61
Le meilleur modèle en prenant en compte the Parsimony principle est le modèle m1 : AR(1)
# Installer le package forecast si nécessaire
# Charger le package forecast
library(forecast)
# Ajuster un modèle ARIMA automatiquement sur les données stationnaires
fit_arima <- auto.arima(diff_cnt)
# Afficher le résumé du modèle ajusté
summary(fit_arima)
## Series: diff_cnt
## ARIMA(0,1,3)
##
## Coefficients:
## ma1 ma2 ma3
## -0.5378 -0.3000 -0.1367
## s.e. 0.0367 0.0414 0.0361
##
## sigma^2 = 1664577: log likelihood = -6161.37
## AIC=12330.74 AICc=12330.8 BIC=12349.05
##
## Training set error measures:
## ME RMSE MAE MPE MAPE MASE
## Training set -44.02235 1286.591 924.9182 97.16658 163.5412 0.8277705
## ACF1
## Training set 0.001903797
residuals_adjusted_deseasonal <- residuals(fit_arima)
plot(residuals_adjusted_deseasonal)
#mean=0 var=cte
par(mfrow = c(2, 1))
acf(residuals_adjusted_deseasonal, main = "ACF of Residuals")
pacf(residuals_adjusted_deseasonal, main = "PACF of Residuals")
#no significant peak on the ACF and PACF
Box.test(residuals_adjusted_deseasonal, type = "Ljung-Box")
##
## Box-Ljung test
##
## data: residuals_adjusted_deseasonal
## X-squared = 0.0026169, df = 1, p-value = 0.9592
aic_value_adjusted<- AIC(fit_arima)
cat("AIC du modèle auto arima:", aic_value_adjusted, "\n")
## AIC du modèle auto arima: 12330.74
Le modèle choisi plus haut, AR(1) donne un AIC légèrement plus élevé, de plus, il est plus simple. et d’après le principe de parsimony, il est préférable de choisir le plus simple.
library(forecast)
forecasts <- forecast(model_m1_deseasonal, h=25)
plot(forecasts)
cat(length(day_data$cnt))
## 731
# Diviser les données en ensembles d'entraînement et de test
train_data <- window(day_data$cnt, start = 1, end = 700)
test_data <- window(day_data$cnt, start = 701)
train_data_ts <- ts(train_data)
test_data_ts <- ts(test_data)
ts.plot(train_data_ts)
diff_train <- diff(train_data_ts, frequency=12)
plot(diff_train)
pacf(diff_train, main="PACF")
acf(diff_train, main="ACF")
On remarque ACF qui décroit donc AR(5) PACF décroît donc MA(3)
auto_model <- auto.arima(train_data)
summary(auto_model)
## Series: train_data
## ARIMA(3,1,1)
##
## Coefficients:
## ar1 ar2 ar3 ma1
## 0.3160 -0.0401 -0.0916 -0.8718
## s.e. 0.0437 0.0410 0.0409 0.0237
##
## sigma^2 = 831740: log likelihood = -5754.54
## AIC=11519.08 AICc=11519.16 BIC=11541.82
##
## Training set error measures:
## ME RMSE MAE MPE MAPE MASE
## Training set 32.40763 908.7347 637.0445 -44.80348 59.12826 0.876924
## ACF1
## Training set -0.006730811
# Ajuster un modèle ARIMA manuellement sur l'ensemble d'entraînement
manual_model <- arima(train_data, order = c(0, 0, 3), seasonal = list(order = c(4, 1, 0), period = 12))
# Ajuster un modèle ARIMA automatiquement sur l'ensemble d'entraînement
auto_model <- auto.arima(train_data)
# Prévoir les prochaines 25 observations avec les deux modèles
manual_forecast <- forecast(manual_model, h = 25)
auto_forecast <- forecast(auto_model, h = 25)
# Plot des prévisions et des données réelles
plot(test_data, col = "blue", ylim = range(test_data, manual_forecast$mean, auto_forecast$mean),
main = "ARIMA Forecast Comparison", xlab = "Time", ylab = "Value")
lines(manual_forecast$mean, col = "red", lty = 2, lwd = 2)
lines(auto_forecast$mean, col = "green", lty = 2, lwd = 2)
# Ajouter une légende
legend("topright", legend = c("Observed", "Manual Forecast", "Auto Forecast"),
col = c("blue", "red", "green"), lty = c(1, 2, 2), lwd = 2)